#include "cells.hpp"
#include "mesh.hpp"
#include "util.hpp"
#include "mesh_testFunctions.hpp"

#include <catch2/catch_test_macros.hpp>
#include <catch2/catch_approx.hpp>
#include <cmath>

using namespace hms;
using Catch::Approx;

TEST_CASE( "Vertices of cells" )
{
	SECTION( "rectangular mesh" )
	{
		auto mesh { testMeshRectangles(true) };
		
		bool ordered {true};
		REQUIRE( exactlyEqual(mesh.cells().vertices( 0, ordered), Array4i{0, 3, 4, 1}) );
		REQUIRE( exactlyEqual(mesh.cells().vertices( 1, ordered), Array4i{1, 4, 5, 2}) );
		REQUIRE( exactlyEqual(mesh.cells().vertices( 2, ordered), Array4i{3, 6, 7, 4}) );
		REQUIRE( exactlyEqual(mesh.cells().vertices( 3, ordered), Array4i{4, 7, 8, 5}) );
		REQUIRE( exactlyEqual(mesh.cells().vertices( 4, ordered), Array2i{6, 7}) );
		REQUIRE( exactlyEqual(mesh.cells().vertices( 5, ordered), Array2i{7, 8}) );
		REQUIRE( exactlyEqual(mesh.cells().vertices( 6, ordered), Array2i{0, 3}) );
		REQUIRE( exactlyEqual(mesh.cells().vertices( 7, ordered), Array2i{3, 6}) );
		REQUIRE( exactlyEqual(mesh.cells().vertices( 8, ordered), Array2i{0, 1}) );
		REQUIRE( exactlyEqual(mesh.cells().vertices( 9, ordered), Array2i{1, 2}) );
		REQUIRE( exactlyEqual(mesh.cells().vertices(10, ordered), Array2i{2, 5}) );
		REQUIRE( exactlyEqual(mesh.cells().vertices(11, ordered), Array2i{5, 8}) );
	}
	SECTION("triangular mesh")
	{
		auto mesh { testMeshTriangles(true) };
		
		bool ordered {true};
		REQUIRE( exactlyEqual(mesh.cells().vertices( 0, ordered), Array3i{0, 5, 4}) );
		REQUIRE( exactlyEqual(mesh.cells().vertices( 1, ordered), Array3i{0, 4, 3}) );
		REQUIRE( exactlyEqual(mesh.cells().vertices( 2, ordered), Array3i{0, 3, 2}) );
		REQUIRE( exactlyEqual(mesh.cells().vertices( 3, ordered), Array3i{1, 0, 2}) );
		REQUIRE( exactlyEqual(mesh.cells().vertices( 4, ordered), Array2i{0, 5}) );
		REQUIRE( exactlyEqual(mesh.cells().vertices( 5, ordered), Array2i{1, 0}) );
		REQUIRE( exactlyEqual(mesh.cells().vertices( 6, ordered), Array2i{1, 2}) );
		REQUIRE( exactlyEqual(mesh.cells().vertices( 7, ordered), Array2i{2, 3}) );
		REQUIRE( exactlyEqual(mesh.cells().vertices( 8, ordered), Array2i{3, 4}) );
		REQUIRE( exactlyEqual(mesh.cells().vertices( 9, ordered), Array2i{5, 4}) );
	}
	SECTION("mixed mesh")
	{
		auto mesh { testMeshMixed(true) };
		
		bool ordered {true};
		REQUIRE( exactlyEqual(mesh.cells().vertices( 0, ordered), Array3i{1, 2, 0}) );
		REQUIRE( exactlyEqual(mesh.cells().vertices( 1, ordered), Array4i{3, 4, 2, 1}) );
		REQUIRE( exactlyEqual(mesh.cells().vertices( 2, ordered), Array2i{2, 0}) );
		REQUIRE( exactlyEqual(mesh.cells().vertices( 3, ordered), Array2i{1, 0}) );
		REQUIRE( exactlyEqual(mesh.cells().vertices( 4, ordered), Array2i{3, 1}) );
		REQUIRE( exactlyEqual(mesh.cells().vertices( 5, ordered), Array2i{4, 2}) );
		REQUIRE( exactlyEqual(mesh.cells().vertices( 6, ordered), Array2i{3, 4}) );
	}
}

TEST_CASE( "Cell edges" )
{
	SECTION( "rectangular mesh" )
	{
		auto mesh { testMeshRectangles(true) };
		
		REQUIRE(sameEntries(mesh.cell(0).edges(), Array4i{0,2,6, 8}));
		REQUIRE(sameEntries(mesh.cell(1).edges(), Array4i{0,3,9,10}));
		REQUIRE(sameEntries(mesh.cell(2).edges(), Array4i{1,2,4, 7}));
		REQUIRE(sameEntries(mesh.cell(3).edges(), Array4i{1,3,5,11}));
		REQUIRE( mesh.cell( 4).edge() ==  4 );
		REQUIRE( mesh.cell( 5).edge() ==  5 );
		REQUIRE( mesh.cell( 6).edge() ==  6 );
		REQUIRE( mesh.cell( 7).edge() ==  7 );
		REQUIRE( mesh.cell( 8).edge() ==  8 );
		REQUIRE( mesh.cell( 9).edge() ==  9 );
		REQUIRE( mesh.cell(10).edge() == 10 );
		REQUIRE( mesh.cell(11).edge() == 11 );
	}
	SECTION("triangular mesh")
	{
		auto mesh { testMeshTriangles(true) };
		
		REQUIRE( sameEntries(
			mesh.cell(0).edges().head(3), Vector3i{1, 3, 8}) );
		REQUIRE( sameEntries(
			mesh.cell(1).edges().head(3), Vector3i{1, 2, 7}) );
		REQUIRE( sameEntries(
			mesh.cell(2).edges().head(3), Vector3i{0, 2, 6}) );
		REQUIRE( sameEntries(
			mesh.cell(3).edges().head(3), Vector3i{0, 4, 5}) );
		REQUIRE( mesh.cell( 4).edge() ==  3 );
		REQUIRE( mesh.cell( 5).edge() ==  4 );
		REQUIRE( mesh.cell( 6).edge() ==  5 );
		REQUIRE( mesh.cell( 7).edge() ==  6 );
		REQUIRE( mesh.cell( 8).edge() ==  7 );
		REQUIRE( mesh.cell( 9).edge() ==  8 );
	}
	SECTION("mixed mesh")
	{
		auto mesh { testMeshMixed(true) };
		
		REQUIRE( sameEntries(
			mesh.cell(0).edges().head(3), Vector3i{0, 1, 2}) );
		REQUIRE( sameEntries(
			mesh.cell(1).edges(), Array4i{0, 3, 4, 5}) );
		REQUIRE( mesh.cell(2).edge() == 1 );
		REQUIRE( mesh.cell(3).edge() == 2 );
		REQUIRE( mesh.cell(4).edge() == 3 );
		REQUIRE( mesh.cell(5).edge() == 4 );
		REQUIRE( mesh.cell(6).edge() == 5 );
	}
}

TEST_CASE( "Cell boundary flags" )
{
	SECTION( "rectangular mesh" )
	{
		auto mesh { testMeshRectangles(true) };
		
		REQUIRE( mesh.cell( 0).isBoundary() == true );
		REQUIRE( mesh.cell( 1).isBoundary() == true );
		REQUIRE( mesh.cell( 2).isBoundary() == true );
		REQUIRE( mesh.cell( 3).isBoundary() == true );
		REQUIRE( mesh.cell( 4).isBoundary() == false);
		REQUIRE( mesh.cell( 5).isBoundary() == false);
		REQUIRE( mesh.cell( 6).isBoundary() == false);
		REQUIRE( mesh.cell( 7).isBoundary() == false);
		REQUIRE( mesh.cell( 8).isBoundary() == false);
		REQUIRE( mesh.cell( 9).isBoundary() == false);
		REQUIRE( mesh.cell(10).isBoundary() == false);
		REQUIRE( mesh.cell(11).isBoundary() == false);
	}
	SECTION("triangular mesh")
	{
		auto mesh { testMeshTriangles(true) };
		
		REQUIRE( mesh.cell( 0).isBoundary() == true );
		REQUIRE( mesh.cell( 1).isBoundary() == true );
		REQUIRE( mesh.cell( 2).isBoundary() == true );
		REQUIRE( mesh.cell( 3).isBoundary() == true );
		REQUIRE( mesh.cell( 4).isBoundary() == false);
		REQUIRE( mesh.cell( 5).isBoundary() == false);
		REQUIRE( mesh.cell( 6).isBoundary() == false);
		REQUIRE( mesh.cell( 7).isBoundary() == false);
		REQUIRE( mesh.cell( 8).isBoundary() == false);
		REQUIRE( mesh.cell( 9).isBoundary() == false);
	}
	SECTION("mixed mesh")
	{
		auto mesh { testMeshMixed(true) };
		
		REQUIRE( mesh.cell(0).isBoundary() == true  );
		REQUIRE( mesh.cell(1).isBoundary() == true  );
		REQUIRE( mesh.cell(2).isBoundary() == false );
		REQUIRE( mesh.cell(3).isBoundary() == false );
		REQUIRE( mesh.cell(4).isBoundary() == false );
		REQUIRE( mesh.cell(5).isBoundary() == false );
		REQUIRE( mesh.cell(6).isBoundary() == false );
	}
}

TEST_CASE( "Ghost cell flags" )
{
	SECTION( "rectangular mesh" )
	{
		auto mesh { testMeshRectangles(true) };
		for ( int i{0}; i<12; ++i)
			REQUIRE( mesh.cell(i).isGhost() == (i >= 4) );
	}
	SECTION( "triangular mesh" )
	{	
		auto mesh { testMeshTriangles(true) };
		for ( int i{0}; i<10; ++i)
			REQUIRE( mesh.cell(i).isGhost() == (i >= 4) );
	}
	SECTION( "mixed mesh" )
	{
		auto mesh { testMeshMixed(true) };
		for ( int i{0}; i<7; ++i)
			REQUIRE( mesh.cell(i).isGhost() == (i >= 2) );
	}
}

TEST_CASE( "Neighbours" )
{
	auto testGhostNeighbours = [](const Mesh& mesh, const auto& neighbours){
		for ( int i{0}; i<mesh.nGhostCells(); ++i )
			REQUIRE( mesh.cell( i+mesh.nDomainCells() ).neighbours()[0] ==
				neighbours(i) );
	};
	
	SECTION( "rectangular mesh" )
	{
		auto mesh { testMeshRectangles(true) };
		
		SECTION( "domain cell neighbours" )
		{
			REQUIRE( sameEntries(
				mesh.cell(0).neighbours(), Array4i{ 1, 2, 6,  8 } ) );
			REQUIRE( sameEntries(
				mesh.cell(1).neighbours(), Array4i{ 0, 3, 9, 10 } ) );
			REQUIRE( sameEntries(
				mesh.cell(2).neighbours(), Array4i{ 0, 3, 4,  7 } ) );
			REQUIRE( sameEntries(
				mesh.cell(3).neighbours(), Array4i{ 1, 2, 5, 11 } ) );
		}
		SECTION( "ghost cell neighbours" )
		{
			ArrayNi<8> neighbours { 2, 3, 0, 2, 0, 1, 1, 3 };
			testGhostNeighbours(mesh, neighbours);
		}
	}
	SECTION( "triangular mesh" )
	{
		auto mesh { testMeshTriangles(true) };
		
		SECTION( "domain cell neighbours" )
		{
			REQUIRE( sameEntries(
				mesh.cell(0).neighbours(), Array3i{ 1, 4, 9 } ) );
			REQUIRE( sameEntries(
				mesh.cell(1).neighbours(), Array3i{ 0, 2, 8 } ) );
			REQUIRE( sameEntries(
				mesh.cell(2).neighbours(), Array3i{ 1, 3, 7 } ) );
			REQUIRE( sameEntries(
				mesh.cell(3).neighbours(), Array3i{ 2, 5, 6 } ) );
		}
		SECTION( "ghost cell neighbours" )
		{
			ArrayNi<6> neighbours {0, 3, 3, 2, 1, 0};
			testGhostNeighbours(mesh, neighbours);
		}
	}
}

TEST_CASE( "Area calculation" )
{
	SECTION( "rectangular mesh" )
	{
		auto mesh { testMeshRectangles(true) };
		
		for ( int i{0}; i<mesh.nCells(); ++i )
			REQUIRE( mesh.cell(i).area() == Approx(1) );
	}
	SECTION( "triangular mesh" )
	{
		auto mesh { testMeshTriangles(true) };
		
		REQUIRE( mesh.cell( 0).area() == Approx( 3) );
		REQUIRE( mesh.cell( 1).area() == Approx( 3) );
		REQUIRE( mesh.cell( 2).area() == Approx( 4) );
		REQUIRE( mesh.cell( 3).area() == Approx( 4) );
		#ifdef MESH_SQUARE_GHOSTS
		REQUIRE( mesh.cell( 4).area() == Approx( 4) );
		REQUIRE( mesh.cell( 5).area() == Approx( 8) );
		REQUIRE( mesh.cell( 6).area() == Approx(10) );
		REQUIRE( mesh.cell( 7).area() == Approx( 8) );
		REQUIRE( mesh.cell( 8).area() == Approx( 4) );
		REQUIRE( mesh.cell( 9).area() == Approx(10) );
		#else
		REQUIRE( mesh.cell(4).area() == Approx(3) );
		REQUIRE( mesh.cell(5).area() == Approx(4) );
		REQUIRE( mesh.cell(6).area() == Approx(4) );
		REQUIRE( mesh.cell(7).area() == Approx(4) );
		REQUIRE( mesh.cell(8).area() == Approx(3) );
		REQUIRE( mesh.cell(9).area() == Approx(3) );
		#endif
	}
	SECTION( "mixed mesh" )
	{
		auto mesh { testMeshMixed(true) };
		
		REQUIRE( mesh.cell( 0).area() == Approx( 2) );
		REQUIRE( mesh.cell( 1).area() == Approx( 4) );
		#ifdef MESH_SQUARE_GHOSTS
		REQUIRE( mesh.cell( 2).area() == Approx( 5) );
		REQUIRE( mesh.cell( 3).area() == Approx( 5) );
		REQUIRE( mesh.cell( 4).area() == Approx( 4) );
		REQUIRE( mesh.cell( 5).area() == Approx( 4) );
		REQUIRE( mesh.cell( 6).area() == Approx( 4) );
		#else
		REQUIRE( mesh.cell(2).area() == Approx(2) );
		REQUIRE( mesh.cell(3).area() == Approx(2) );
		REQUIRE( mesh.cell(4).area() == Approx(4) );
		REQUIRE( mesh.cell(5).area() == Approx(4) );
		REQUIRE( mesh.cell(6).area() == Approx(4) );
		#endif
	}
}

TEST_CASE( "Cells::circumf(): Circumference calculation" ){
	SECTION( "rectangular mesh" )
	{
		auto mesh { testMeshRectangles(true) };
		for (int i{0}; i<mesh.nDomainCells(); ++i)
			REQUIRE( mesh.cell(i).circumf() == Approx(4) );
	}
	SECTION( "triangular mesh" )
	{
		auto mesh { testMeshTriangles(true) };
		
		REQUIRE( mesh.cell( 0).circumf() == Approx( 2. + 2*sqrt(10) ) );
		REQUIRE( mesh.cell( 1).circumf() == Approx( 2. + 2*sqrt(10) ) );
		REQUIRE( mesh.cell( 2).circumf() == Approx( 2*(sqrt(10)+sqrt(2) ) ) );
		REQUIRE( mesh.cell( 3).circumf() == Approx( 2*(sqrt(10)+sqrt(2) ) ) );
	}
	SECTION( "mixed mesh" )
	{
		auto mesh { testMeshMixed(true) };
		
		REQUIRE( mesh.cell( 0).circumf() == Approx( 2. + 2*sqrt(5) ) );
		REQUIRE( mesh.cell( 1).circumf() == Approx( 8) );
	}
}

TEST_CASE( "Centroids" )
{
		SECTION( "rectangular mesh" )
	{
		auto mesh { testMeshRectangles(true) };
		SECTION( "Domain cells" )
		{
			REQUIRE( approxEqual( mesh.cell( 0).centroid(), Vector2s{-0.5,-0.5} ) );
			REQUIRE( approxEqual( mesh.cell( 1).centroid(), Vector2s{-0.5, 0.5} ) );
			REQUIRE( approxEqual( mesh.cell( 2).centroid(), Vector2s{ 0.5,-0.5} ) );
			REQUIRE( approxEqual( mesh.cell( 3).centroid(), Vector2s{ 0.5, 0.5} ) );
		}
		SECTION( "Ghost cells" )
		{
			REQUIRE( approxEqual( mesh.cell( 4).centroid(), Vector2s{ 1.5,-0.5} ) );
			REQUIRE( approxEqual( mesh.cell( 5).centroid(), Vector2s{ 1.5, 0.5} ) );
			REQUIRE( approxEqual( mesh.cell( 6).centroid(), Vector2s{-0.5,-1.5} ) );
			REQUIRE( approxEqual( mesh.cell( 7).centroid(), Vector2s{ 0.5,-1.5} ) );
			REQUIRE( approxEqual( mesh.cell( 8).centroid(), Vector2s{-1.5,-0.5} ) );
			REQUIRE( approxEqual( mesh.cell( 9).centroid(), Vector2s{-1.5, 0.5} ) );
			REQUIRE( approxEqual( mesh.cell(10).centroid(), Vector2s{-0.5, 1.5} ) );
			REQUIRE( approxEqual( mesh.cell(11).centroid(), Vector2s{ 0.5, 1.5} ) );
		}
	}
	SECTION( "triangular mesh" )
	{
		auto mesh { testMeshTriangles(true) };
		
		SECTION( "Domain cells" )
		{
			REQUIRE( approxEqual( mesh.cell(0).centroid(), Vector2s{ 1,  1} ) );
			REQUIRE( approxEqual( mesh.cell(1).centroid(), Vector2s{ 0,  2} ) );
			REQUIRE( approxEqual( mesh.cell(2).centroid(), Vector2s{-4.0/3, 4.0/3} ) );
			REQUIRE( approxEqual( mesh.cell(3).centroid(), Vector2s{-5.0/3,-1.0/3} ) );
		}
		SECTION( "Ghost cells" )
		{
			#ifdef MESH_SQUARE_GHOSTS
			REQUIRE( approxEqual( mesh.cell(4).centroid(), Vector2s{ 1, -1} ) );
			REQUIRE( approxEqual( mesh.cell(5).centroid(), Vector2s{ 0, -2} ) );
			REQUIRE( approxEqual( mesh.cell(6).centroid(), Vector2s{-4, -1} ) );
			REQUIRE( approxEqual( mesh.cell(7).centroid(), Vector2s{-3,  3} ) );
			REQUIRE( approxEqual( mesh.cell(8).centroid(), Vector2s{ 0,  4} ) );
			REQUIRE( approxEqual( mesh.cell(9).centroid(), Vector2s{ 3,  2} ) );
			#else
			REQUIRE( approxEqual( mesh.cell(4).centroid(), Vector2s{ 1.0  ,-1.0  } ) );
			REQUIRE( approxEqual( mesh.cell(5).centroid(), Vector2s{-1.0/3,-5.0/3} ) );
			REQUIRE( approxEqual( mesh.cell(6).centroid(), Vector2s{-9.8/3,-2.6/3} ) );
			REQUIRE( approxEqual( mesh.cell(7).centroid(), Vector2s{-8.0/3, 8.0/3} ) );
			REQUIRE( approxEqual( mesh.cell(8).centroid(), Vector2s{ 0.0  , 4.0  } ) );
			REQUIRE( approxEqual( mesh.cell(9).centroid(), Vector2s{ 2.2  , 1.4  } ) );
			#endif
		}
	}
	SECTION( "mixed mesh" )
	{
		auto mesh { testMeshMixed(true) };
		
		SECTION( "Domain cells" )
		{
			REQUIRE( approxEqual( mesh.cell(0).centroid(), Vector2s{ 1,  8./3} ) );
			REQUIRE( approxEqual( mesh.cell(1).centroid(), Vector2s{ 1,  1} ) );
		}
		SECTION( "Ghost cells" )
		{
			#ifdef MESH_SQUARE_GHOSTS
			REQUIRE( approxEqual( mesh.cell(2).centroid(), Vector2s{ 2.5, 3.5} ) );
			REQUIRE( approxEqual( mesh.cell(3).centroid(), Vector2s{-0.5, 3.5} ) );
			REQUIRE( approxEqual( mesh.cell(4).centroid(), Vector2s{-1  , 1  } ) );
			REQUIRE( approxEqual( mesh.cell(5).centroid(), Vector2s{ 3  , 1  } ) );
			REQUIRE( approxEqual( mesh.cell(6).centroid(), Vector2s{ 1  ,-1  } ) );
			#else
			REQUIRE( approxEqual( mesh.cell(2).centroid(), Vector2s{ 6.2/3, 3.2} ) );
			REQUIRE( approxEqual( mesh.cell(3).centroid(), Vector2s{-0.2/3, 3.2} ) );
			REQUIRE( approxEqual( mesh.cell(4).centroid(), Vector2s{-1  , 1  } ) );
			REQUIRE( approxEqual( mesh.cell(5).centroid(), Vector2s{ 3  , 1  } ) );
			REQUIRE( approxEqual( mesh.cell(6).centroid(), Vector2s{ 1  ,-1  } ) );
			#endif
		}
	}
}



#define TEST_CELL_PROPERTIES \
	for ( Index i{0}; i<mesh.nCells(); ++i ){ \
		using hms::approxEqual; \
		using hms::exactlyEqual; \
		REQUIRE( mesh.cells().area   (i) == Approx(mesh.cell(i).area   () ) ); \
		REQUIRE( mesh.cells().circumf(i) == Approx(mesh.cell(i).circumf() ) ); \
		REQUIRE( approxEqual(mesh.cells().centroid(i), mesh.cell(i).centroid() ) ); \
		REQUIRE( mesh.cells().nEdges    (i) == mesh.cell(i).nEdges    () ); \
		REQUIRE( mesh.cells().isGhost   (i) == mesh.cell(i).isGhost   () ); \
		REQUIRE( mesh.cells().isBoundary(i) == mesh.cell(i).isBoundary() ); \
		REQUIRE( exactlyEqual(mesh.cells().neighbours(i), mesh.cell(i).neighbours() ) ); \
		REQUIRE( exactlyEqual(mesh.cells().vertices(i), mesh.cell(i).vertices() ) ); \
		if ( i<mesh.nDomainCells() ){ \
			REQUIRE( exactlyEqual(mesh.cells().edges(i), mesh.cell(i).edges() ) ); \
		} else { \
			REQUIRE( mesh.cells().edge (i) == mesh.cell(i).edge () ); \
		} \
	}


TEST_CASE( "Wrapper object Cell" ){
	SECTION( "rectangular mesh" ){
		auto mesh { testMeshRectangles(true) };
		TEST_CELL_PROPERTIES
	}
	SECTION( "triangular mesh" ){
		auto mesh { testMeshTriangles(true) };
		TEST_CELL_PROPERTIES
	}
	SECTION( "mixed mesh" ){
		auto mesh { testMeshMixed(true) };
		TEST_CELL_PROPERTIES
	}
}